Toroidal L and H equilibria with axisymmetric rotations 
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ABSTRACT 



Axisymmetric toroidal equilibria with toroidal and poloidal rotations are 
solved with a specific set of source functions. The two independent solutions 
are associated to L and H modes. The L/H transition is regarded as a bifurca- 
tion from one equilibrium configuration to another, under strong external heating 
and pellet injection to shape temperature and density profiles. Because of the 
steep edge gradient of the H solution, large static radial electric field, zonal fiow, 
and improved confinements, come as consequences, not causes, of the H mode. 

Subject headings: Toroidal Equilibria 
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I. INTRODUCTION 

The discovery of the H mode with enhanced confinement in ASDEX [1] had openned 
a new age of tokamak fusion research. Extensive experimental works were performed to 
identify the key procedures and signatures that carried the L to the H confinement mode 
[2,3]. These procedures included pellet injection for density enhancement, neutral beam or 
radio frequency heating for temperature and electric conductivity profiling, etc [4,5]. The 
experimental signatures of L/H transition were a drastic reduction of the hydrogen 
emission and a sudden decrease of the plasma fioating potential. These were accompanied 
by a good enhancement of plasma density and energy confinements, plus a steepening of 
the edge plasma profile and a reduction of magnetohydrodynamic (MHD) activities [6,7]. 
Observationally, the L/H transition was accompanied by an increase in the toroidal and 
poloidal rotations [3]. Nevertheless, it was not clear that this velocity increase was the 
cause or was the consequence of the transition. 

The L/H transition is believed to be caused by divertor shaping of the edge plasma, 
under high levels of heating power. Action of the divertor plus the edge gradient pump 
up the static radial electric field that drives a zonal fiow near the edge rational magnetic 
surface [8,9]. This zonal flow is thought bo cause a transport barrier that enhances density 
and energy confinements. The plasma is thought to self-organize gradually under this 
scenario to reach the H mode [10]. Here, we take a different approach to view the L/H 
transition. We consider the toroidal and poloidal rotations on the equilibrium scaling, not 
transport scaling, as equilibrium parameters. We solve for rotational equilibria, under 
specific source functions, in spherical coordinates by seeking toroidal solutions [11]. There 
are two independent solutions for rotational equilibrium, and we associate them to the 
L and H modes. Under this equilibrium configuration approach, the transition is seen 
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as a bifurcation from one equilibrium to another, under external drivings such as pellet 
injection for plasma density and strong heating for temperature profiling. Because of the 
steep edge gradient of the H mode, the large static radial electric field, the zonal flow, and 
the associated confinements, come as the natural consequences, not causes, of the H mode 
solution. 
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II. ROTATIONAL EQUILIBRIA 

To discuss rotational equilibria, we begin with the time-dependent MHD equations 



-^ + V-{pv) = 0, (1) 

p{^ + {v-W)v} ^ JxB-Wp , (2) 
dB ^ 

— ^ -VxE ^Vx(vxB) , (3) 

VxB = ijj , (4) 

V • S = , (5) 

P = pv] ■ (6) 



Here, p is the mass density, p is the plasma pressure, v is the bulk velocity, J is the current 

— * 

density, B is the magnetic field, Vg is the ion acoustic speed, ji is the free space permeability. 
With axisymmetry, the magnetic field and the current density can be represented by two 
scalar functions in standard spherical coordinates 



A 1 f)P BP 

B = A,{VP X V0 + gV0) = ^{+-^, --f. +Q} , (7) 

rsmfc' r oO or 

^ rsin^'^^ra^' dr ' dr^ r'' de^siuO dO^^ ' ^ ^ 

Here, Aq carries the physical dimension of poloidal magnetic fiux such that P is a 
dimensionless function. Also, we can write the axisymmetric poloidal and toroidal rotations 

as 



(9) 
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Likewise, carries the physical dimension of poloidal velocity flux such that P' is a 
dimensionless function. With axisymmetry and incompressible fluid condition, V • w = 0, 
steady state in Eq. [1] requires 



Vp-v = A^Vp ■ (VP' X V0) = A[,(Vp X VP') ■ V0 = , 

p = PoPin , (10) 

which requires the poloidal velocity and the mass density have the same level contours. 
Here, po carries the dimension and amplitude of mass density, and p(P') is a dimensionless 
function. As for Eq. [3l by Eq. [7] and Eq. [9], we note that v x B would be null and steady 
state in Eq. [3] would be warrented with 



P' = aP , (11a) 
Q' = aQ . (lib) 



As a result, the velocity field and the magnetic field are parallel, generating an emf-free 
velocity field 



V = a^B = gB . (12) 

^0 



With Eq. m and Eq. [6], toroidal plasma equilibria of Eq. [2] with full axisymmetric rotations 
are described by 
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(1 - iipog''p(P))(V xB)xB- iipovyp(P) = ^i2pog^p(P)VB^ = . (13) 

We consider the rotational scalar pressure be much smaller than the plasma pressure, 
thereby giving the second equality in the above equation. 
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III. ROTATIONAL GRAD-SHAFRANOV EQUATION 

We seek to solve Eq. [13] for toroidal solutions. This equation renders three components. 
By axisymmetry, the component contains only the magnetic force, and it is 

dPdQ_dPdQ^ 

dr dd de dr ' ^ ' 

Q{T,e) = Q{P{r,e)) . (14b) 



As for the 6 component, it reads 



+ -) + -^} = +(-)V=m^9-l„(l-Wp(i')) • (15) 

This equation is the rotational counterpart of the Grad-Shafranov equation of axisymmetric 
toroidal plasma equilibrium, represented in spherical coordinate system. The three 
terms on the left side represent the nonlinear force-free field with fiJ = K{P)B, where 
K{P) = dQ/dP is a scalar function. This can be verified from Eq. [7] and Eq. [8] when we 
impose /iJ = K{P)B. In particular, we would have the linear force- free field should we take 
Q^{P) = {dPy with constant K{P) = a. The term on the right side is the plasma pressure 
balance. The magnetic function Q^{P) and the mass density p{P) are source functions that 
need to be specified. This second order partial differential equation has two independent 
solutions. Finally, the r component of Eq. [13] reads 



d 



-2 ^^^^ de^sine de 

d 



2 dP 



+ (^)Vsin2 0— (1- 



(16) 
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Comparing Eq. [16] to Eq. [151 we note that these two equations are identical. The r 
component is simply the self-consistent condition of the 9 component. 

To solve Eq. [TB] analytically, we take the source functions as 



(17a) 
(17b) 



1 - e 



-pi 



Writing P{r,6) = R{r)Q{6), the rotational Grad-Shafranov equation reads 



= -^qiROy-^rUin^e = -^^q(Rey-^r^ sin^ 6 . 
The variables of this equation could be separated by taking g = 1 to give 



>(PR 



+ [{arf - n{n + 1)]R 



(v4Qa' 



2\2 



e 



-Ai(ar) 



e 



, (19a) 
, (19b) 



where we have denoted x = cos6, and used n{n + 1) as the separation constant. The factor 
v'^ / {AqO^)"^ is proportional to v'^/v'^^i, where Vpoi is the average poloidal rotational velocity. 
The first equation gives 
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e(x) = (l-x^)^^l^ = (l-x^) , (20) 
dx 

where Pn{x) is the Legendre polynomial. We have taken n = 1 to get the second equality. 
As for the second equation, the 9 dependent part on the right side disappears by having 
n = 1. The solution is given by R{r) = Ro{r) + Ri{r), where Ro{r) and Ri{r) are the 
homogeneous and particular solutions. The homogeneous solution is described by 



Ro{r) = arjn{ar) + \oaryn{ar) , (21) 

where jn{z) and yn{z) are the oscillating spherical Bessel functions, and Aq is a constant. 
Together with Aq defined in Eq. [71 there are two constants for Roir). As for the particular 
solution, we have 



z^^ + [z^-n{n + l)]R, = -A.z" , 

R^[r) = -Aiiarf = -Aiz^ . (22) 

We note that the homogeneous solutions, Ro{r) and Q{x), correspond to the linear or 
nonlinear force-free solutions of the left side of Eq. [151 The plasma pressure term on 
the right side appears only in the particular solution, Ri{r), that keeps the pressure 
balance. The homogeneous radial solution is an oscillating function in z = ar, which has 
sucessive maxima, and the homogeneous meridian solution has a lobe peaked at x = 0. 
The superposition of the particular radial solution only slightly modifies the homogeneous 
solutions. We could use the region between z = and the first root of jn{z), with n = 1, to 
describe low aspect ratio high (3 toroidal plasma equilibria. 
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IV. MAGNETIC AND CURRENT STRUCTURES 

With the spatial structure solved, the magnetic field components are given by 



Br 



IdP 



Bn 



r sin 9 r 89 
1 dP _ ldR{r) 

r sin 9 Or r dr 

Brh 



r 



1 



.2,2)1/2 
1 

r sin 9 



dx 
Q(x) 



Q{P) 



(23a) 
(23b) 
(23c) 



The solution R{r) vanishes at some r where we have Br{r) = 0. The solution 0(x) also 
vanishes at some x. Together they describe the magnetic fields. Within this region of (r, x), 
the topological center defined by dR{r)/dr = and dQ{x)/dx = has Br = and Bq = 0. 
This is the magnetic axis, r = r^, where the magnetic field is entirely toroidal. The field 
lines about this center are given by 



B^ ^ B^ ^ B^ 
dr rd9 r sin f-^-^ 



By axisymmetry, the third group is decoupled from the first two groups. For the field lines 
on an (r — 9) plane, we consider the first equality between Br and Bg which gives 



P(r,x) = R{r)e{x) = C . (25) 

The nested poloidal field lines are given by the contours of P{r,x) on the (r — x) plane. 
At the topological center, we have Q{x) maximum and R{r) maximum, so that P{r,x) is 
maximum. Since r sin 9 is the distance of a point on the (r — x) plane to the z axis, Eq. I23cl 
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states that the hne integral of around the circle on the azimuthal plane is measured by 
2nQ, 



27rr sin 6B^ = 2'kQ = jjl^ . 

This line integral about the axis of symmetry is maximum at the topological center. Also, 
it is evident that Q is equivalent to the axial current, where the constant part Qo amounts 
to a uniform component. As for P, we evaluate the poloidal magnetic flux by integrating 
Eq. I23bl on the x = plane over a cross section to give 




27irBedr = -2n{P{z) - P{z^)) . 



As for the current density of Eq. [8], making use of the Grad-Shafranov equation of Eq. [T5] 
gives 

l^J = ^A+-% +ia'P - ^i^fr' sin^ 6^ ln(l - Wp(P))} • (26) 

r sm fc' r o6 or Aq oP 

Analogous to the magnetic field lines, the current density field lines are given by 

Jr_ _ Je _ J(f> /2^\ 
dr rdO r sin Odcf) 

Considering the first equality, the poloidal current density contours are given by 



Q{r,x) = C 



(28) 
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V. L AND H MODES 

We note that there are two independent solutions for Ro{r) in Eq. [211 The first one is 
zji{z) which vanishes at z = 0. With zi and Z2 as the first and second zeros, the region 
bounded by < z < could be used to describe spheromak and high (3 low aspect ratio 
tokamak equilibria. The second one is zyi{z) which diverges at z = 0. Since our domain of 
interest in tokamak plasmas excludes z = 0, the singularity of yi is irrelevant. The region 
bounded hj zi < z < Z2 could also be used to describe high (3 low aspect ratio tokamak 
equilibria. The functions zji{z) and zyi{z) are shown in Fig.l. The poloidal magnetic 
contours of Eq. [25] for zji{z) in the interval < z < Zi are shown in Fig.2. In particular, 
this solution could also be applied to spheromaks where 2; = is accessible to plasma 
equilibria. Similar contours for zi < z < Z2, Z2 < z < Z3, and so on, can be obtained 
to represent tokamak plasmas of different aspect ratios. The contours for zyi{z) in the 
interval zi < z < Z2 are shown in Fig. 3. In order to illustrate the essential features, we have 
neglected the particular solution -Ri(r), and have taken R{r) = Ro{r). The contour levels 
are taken at 0.95,0.9,0.7,0.5,0.3,0.1 of the respective peak value. The external contours 
indicate high poloidal fields, and internal contours for low poloidal fields. These contours 
also indicate the poloidal rotations with rotation velocity high on the outside and low on 
the inside. 

Bounded by a smaller interval zi < z < Z2, we note that Fig. 3 of zyi{z) has a more 
localized domain and steeper edge profile than the equilibrium in Fig.2, which is described 
by zji{z) in the larger interval < z < zi. Including the negative valued particular 
solution -Ri(r) further steepens the edge gradient. We associate zji{z) and zyi{z) to the 
L and H mode respectively. We also note that the poloidal and toroidal magnetic fields 
are plotted in normalized radial coordinate z = ar. In laboratory plasmas, the fields are 
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measured in terms of radius r. To connect to our normalized results, we need to determine 
the normalizing parameter a. This can be done by considering the magnetic axis of a 
laboratory plasma, say in the zji{z) mode bounded by < z < Zi, through a^r* = Zi, = 2.7. 
Defined by the divertor scrape-off, the radial range of laboratory plasma, Va < r < r^, can 
now be converted to < Za < z < Zb < zi with zi = 4.5. In the case of zyi{z) mode, 
we have Qyr^ = z^ = 4.5, giving ay/aj = 4.5/2.7. The radial range can be converted to 
Zi < Za < z < Zb < Z2 with zi = 2.8 and Z2 = 6.1. Experimentally, this change of the 
normalizing parameter from aj to ay could be accomplished by pellet injection and high 
power external heating. 

In the L mode, the zji{z) profile is more diffused between a larger interval of the 
zeros. The divertor action removes the edge plasma to the Za < z < Zb domain, with 
poloidal rotation velocity contours corresponding to such domain. In the H mode, due to 
the compactness of the interval between zeros of zyi{z) profile, the plasma equilibrium fits 
within the toroidal machine vessel naturally with much less divertor shaping. The plasma 
equilibrium occupies probably the entire domain zi < z < Z2, or a large part of it. As a 
result, poloidal rotation contours of the H mode cover not just the central part but also the 
high velocity part on the outside. By going from L to H mode, the rotation contours within 
the plasma cross-section, defined by the divertor action, are enlarged from a partial central 
profile to an almost complete profile. Observed at a fixed position at the plasma edge, we 
would have the impression that the rotation velocity has been speeded up. 

The corresponding mass density contours of Eq.(17b) for L and H modes are shown in 
Fig. 4 and Fig.5 respectively. The profiles along x = horizontal cut are shown in Fig. 6. 
Although the two modes are presented in one same figure showing approximately the same 
dimensionless amplitudes, the physical amplitude and dimension is given by po defined in 
Eq. [TOl As a result, the mass density of the H mode could be much larger than that of the 
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L mode. The essence of Fig. 6 is to show the relative shape of the mass density profiles for 
the two modes. We have suggested the identification of L mode to the zji (z) solution in the 
{0,Zi) domain, and H mode to the zyi{z) solution in the (^1,^2) domain. To map these z 
domains to the same r domain of machine vessel, We have used two different normalizing 
parameters aj and ay for the source function Q'^{P) of Eq. I17al Since ay = 1.7a this would 
require a substantial increase of toroidal magnetic field according to Eq. I23c[ To avoid this 
substantial toroidal field enhancement, the H mode could be generated by superimposing 
the zyi{z) solution to the zji^z) solution, without displacing significantly the z domain. 
As an example, with Aq = 1, the profile of {zji{z) + zyi{z)) is shown in Fig. 7 indicating 
zi = 1.9 and Z2 = 5.3 with a maximum at = 3.6. 
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VI. DISCUSSIONS AND CONCLUSIONS 

Experimentally, due to the divertor action on the edge plasma, there is an electrostatic 
field normal to the magnetic surfaces. This field is particularly large in the H mode 
configuration because of the steep edge gradient. Since the poloidal field lines are described 
by the P contours, this field can be written as ^ = — V$(P), which warrants dB/dt — 
for equilibrium. Interacting with magnetic islands on a rational surface, this electric field 
drives zonal fiows that establish transport barriers for better plasma density and energy 
confinements. With the rotational toroidal equilibrium approach, the L/H transition 
amounts to a bifurcation of the L equilibrium to the H equilibrium, under the actions of 
external pumping through pellet injection and strong auxiliary heating for density and 
temperature profile shaping. The normal electric field, zonal fiows, and transport barriers, 
come as consequences, not cuases, of the steep edge gradient of the H mode. 

We have solved toroidal plasma equilibria with axisymmetric toroidal and poloidal 
rotations, that are self-similar to the corresponding magnetic fields. The rotational 
Grad-Shafranov equation in spherical coordinates is solved for toroidal solutions, under the 
assumption that the scalar rotational pressure is much less than the plasma pressure. With 
a specific set of source functions, there are two independent homogeneous radial modes 
given by zji{z) and zyi{z). The zji{z) mode in the region < z < zi and the zyi{z) mode 
in the region zi < z < Z2 could be applied to current large scale tokamaks. The zji{z) 
mode has a diffuse edge profile over a larger z domain, and the zyi{z) mode has a steep 
edge profile over a smaller z domain. We associate them to the L and H modes respectively. 
The L/H transition amounts to a bifurcation of one equilibrium configuration to another, 
following a change of the normalizing parameter from aj to ay. Experimentally, this change 
of parameter could be achived by pellet injection and large external heating. 
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Fig. 1. — The functions zji{z), in thin hne, and zyi{z), in thick hne, of the two independent 
homogeneous Rq{z) radial solutions are plotted as a function of z. 
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Fig. 2. — The poloidal magnetic field lines of Ro{z) = zji{z) L mode are shown in polar 
plots with increasing contour labels from outer to inner contours. These contours are shared 
by the poloidal rotation. The axes are labelled in z together with x = cos 9. 
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Fig. 3. — The poloidal magnetic field lines of Ro{z) = zyi{z) H mode are shown in polar 
plots with increasing contour labels from outer to inner contours. 




Fig. 4. — The mass density contours of Ro{z) = zji{z) L mode are shown in polar plots with 
increasing contour labels from outer to inner contours. 




Fig. 5. — The mass density contours of Ro{z) = zyi{z) R mode are shown in polar plots 
with increasing contour labels from outer to inner contours. 
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Fig. 6. — The mass densities along a; = for zji{z) mode, in thin hne, and zyi{z) mode, in 
thick hne, are plotted as a function of z. 




Fig. 7. — The function ,2^1(2;) + zyi{z) is plotted as a function of z indicating zi = 1.9 and 
Z2 = 5.3. 



